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ABSTRACT 

The structure of static atmospheres around unmagnetized neutron stars 
undergoing steady, spherical accretion is discussed. We focus on the “hot” 
configurations presented by Turolla et al. (1994) and calculate the radiation 
spectrum using a characteristics method. In particular, it is found that e + -e~ 
pair production may affect significantly the external atmospheric layers, where 
positron and proton number densities become of the same order. The consequent 
increase of the scattering opacity lowers the Eddington limit and this, in turn, 
may drive a dynamical instability if the accretion luminosity is large enough, 
ultimately producing a rapid expulsion of the envelope. If “hot” states are 
indeed accessible, this mechanism could give rise to transient phenomena in 
hard X-rays of potential great astrophysical interest. 

Subject headings: accretion, accretion disks — radiative transfer — stars: 
neutron — X-rays: general 



1. Introduction 


The problem of calculating the radiation spectrum produced in the atmosphere of an 
unmagnetized, accreting neutron star (NS) received wide attention in the past, starting 
from the pioneering work by Zel’dovich, & Shakura (1969, hereafter ZS) where the case 
of stationary, spherical accretion was first discussed. In the commonly accepted picture, 
steaming from ZS’ original analysis and corroborated by the more detailed numerical 
work by Alme, & Wilson (1973), the emerging spectrum is essentially a blackbody plus 
a high-energy Compton tail. It can be characterized by two parameters, the accretion 
luminosity (or, equivalently, the accretion rate) and the column density yo corresponding 
to the penetration range of incoming protons in the NS atmosphere. Recently Zarnpieri 
et al. (1995, paper II) extended ZS’s results to cover the low luminosities typical of old, 
isolated NSs accreting the interstellar medium (L <; 10 34 erg/s), showing that deviations 
from a blackbody become more conspicuous as the luminosity decreases. 

Turolla et al. (1994, paper I) re-examined the issue of the thermal and radiative 
properties of static atmospheres around accreting NSs and pointed out that a new class of 
much hotter equilibrium solutions may exist if the accretion luminosity exceeds a certain 
critical value, in analogy with what has been already found in black hole accretion (see e.g. 
Park 1990; Nobili, Turolla, & Zarnpieri 1991). The existence of two states is related to the 
different role thermal and non-thermal processes may play in exchanging energy between 
photons and electrons in the atmosphere: in the “cold” state (a la ZS) the energy released 
by accretion is radiated away essentially via bremsstrahlung, whereas in the “hot” state 
Compton scattering dominates (in this respect see also Burger, & Katz 1980). 

In paper I “hot” solutions were found solving the frequency-integrated transfer 
problem. The lack of any information about the angular and frequency dependence of the 
radiation field made impossible to evaluate the photon-photon pair production, which was 
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nevertheless recognized to be potentially important, and led to introduce rather drastic 
assumptions in the calculation of the Compton energy exchange rate. 

In this paper we try to overcome some of these limitations, in the attempt to shed light 
on the relevance of “hot” solutions. To get some physical insight without resorting to a 
full time-dependent calculation, we consider the stationary case and analyze the effects of 
different processes separately. First the radiation field is calculated using the characteristics 
method introduced by Zane et al. (1996, paper III). Pair production is neglected and the 
frequency-dependent transfer equation is solved coupled to the hydrostatic and energy 
balance. The resulting intensity is then used to estimate the pair number density n + 
produced by photon-photon interactions. Although not fully self-consistent, this approach 
shows that a large pair production is expected in the external atmospheric layers. This 
opens the possibility that such configurations are unstable with respect to the onset of 
a radiatively-driven wind, produced by the increased opacity in the pair-rich plasma. 
Astrophysical applications of “hot” atmospheres are suggested in the context of hard X-ray 
transient phenomena. 


2. The Model 

In this section we briefly outline the input physics of our model, referring to previous 
works (ZS; Alme, & Wilson 1973; paper I and II) for all details, and introduce our radiative 
transfer calculations. 

We consider an unmagnetized, non-rotating neutron star with mass M* = 1M 0 
and radius i?* = 10 6 cm, surrounded by a static atmosphere where the accreting flow 
is decelerated by Coulomb collisions and/or plasma interactions. As the flow hits the 
outermost stellar layers, the proton bulk kinetic energy is converted into electron thermal 
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energy and then re-emitted as electromagnetic radiation. The heat injected per unit time 
and mass, Wh, is given by (see Bildsten, Salpeter, & Wasserman 1992; paper II) 


W h 


Loo _ 1 + Vth/ V i _ 

SnRlyoVg [1 - (1 - vjjvf) (y/y 0 )} 1/2 


( 1 ) 


where y g = (1 — R g /R*) 1 ^ 2 is the gravitational redshift factor, R g = 2GM*/c 2 , 
vl h = 3 KT/rrip, T is the gas temperature, vf = c 2 (1 — L^/Le) ( R g /R *), y = jp pdR is 
the column density and Le is the Eddington luminosity. As discussed in paper I, these 
solutions show no significant envelope expansion, so the radial coordinate can be safely 
assumed equal to the star radius. The atmospheric structure can be solved in plane-parallel 
geometry, using as independent variable the scattering optical depth r = n es y, where n es 
is the Thomson opacity. In this case, the luminosity profile can be immediately derived 
from the zero-th frequency-integrated moment equation (see paper II). The run of pressure 
P and temperature T are given by the solution of the hydrostatic and energy balance 
equations (see again paper I and II) 


dP _ GM * / K\ L \ 
dr y 2 g RlK es \ Va K es L E ) 

— = k p faT 4 - ^-u) + (T - A)c , (3) 

c V up J 

with the boundary condition of vanishing pressure at the top of the atmosphere. In equation 
(3) U is the radiation energy density, up, kq and «q are the Planck, absorption, and flux 
mean opacities and (T — A)c is the Compton energy exchange rate. The effective heating 
W represents the fraction of IT),, which is converted into electromagnetic radiation within 
the atmosphere; W coincides with IT/, if no energy is injected at the inner boundary, that is 
if L in = 0. Since at large depths the radiation held becomes more and more isotropic, L tn is 
indeed small, although in a realistic situation it is not exactly zero. On the other hand, the 
thermal structure of both the “cold” and “hot” solutions presented in paper I proved not 
sensitive to the assumed value of the luminosity at the bottom of the atmosphere, provided 
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it is not too large compared with L^. In particular, “hot” atmospheres exhibit a rather 
sharp drop in the gas temperature which separates two distinct regions characterized by 
different thermal properties. In the external layers, where the bulk of the luminosity is 
produced, the energy released by accretion is radiated away via Compton cooling, while in 
the deeper, denser layers LTE is attained and the temperature profile is mainly determined 
by bremsstrahlung equilibrium. For these solutions L in also coincides with the value of the 
radiative flux at the top of the cold region. 

There are two important points concerning the “hot” models that were not thoroughly 
discussed in paper I and deserves further analysis. Both of them concern the possibility of 
establishing and maintaining a rather steep temperature gradient between the hot to the 
cold part of the atmosphere, so that a “hot” solution can indeed set in. First of all we note 
that the position of the temperature drop is practically unrelated the value of the proton 
stopping length. Although for y 0 = 20 (the value we will use here) they practically coincide, 
we stress that, in general, the value of the scattering depth at which the transition occurs 
depends only on the relative efficiency of Compton and bremsstrahlung heating-cooling 
and, as numerical tests show, it is always close to r tr = 10 ( y tr = 25). The presence of a 
sudden decrease in T can be understood comparing the free-free and Compton thermal 
times 

tc_ = € ff (1 - U/aT 4 ) 

t ff Ak\T — T^KeJJ/m e c 2 ’ {> 

where e/f is the free-free emission and T 7 is the radiation temperature. As the scattering 
depth increases T and T 7 become closer, because non-conservative scatterings tend to 
establish thermal equilibrium, and U either increases or stays constant. In the hot, 
effectively thin layers U <C aT 4 and at large enough r, where T ~ T 7 , the time scale for 
free-free cooling becomes shorter than the Compton time. The plasma temperature drops 
until aT 4 ~ U in order for the energy balance to be satisfied and the radiation temperature 
follows. 
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The second important point concerns the efficiency of electron conduction at the 
interface between the hot and cold regions. This effect was not considered in paper I. In 
order to get a quantitative estimate of the relevance of thermoconduction we compare the 
energy flux due to conduction, F c = u c pdT/dy, with the radiative flux, F r = L/AtxK 2 . Using 
the expression for v c given in Zel’dovich, & Raizer (1967) it is 



2.5 x 10" 3 


/ J' \ 5/2 

llOV 



Utr 

Ay 


( 5 ) 


where Ay is the width of the transition layer and T is the temperature in the hot region. 
As the previous equation shows, electron conduction is indeed efficient in transferring heat 
from the hot to the cold region. However, if a fraction of the luminosity is produced below 
y tr , Lin/Loo ~ 0.2 — 0.3 for L <*, ~ 0.1 L E , conduction will limit the temperature of the hot 
region to ~ 10 9 K, smear the jump over Ay ~ 2 — 3 and then cease to be important. The 
stationary temperature profile will then look quite similar to that discussed by Zel’dovich, 
& Raizer (1967) in the context of shock waves in a plasma. Although the real profile will 
deviate from that of paper I (mainly in the absence of the peak at T ~ 10 10 K, see curve b 
in figure 2, which will be leveled if conduction is taken into account), the emerging picture 
seems to be substantially unaltered by the inclusion of thermal conduction. This conclusion 
is reinforced by the result of the frequency-dependent calculation presented here which 
shows that the temperature profile is nearly isothermal at T ~ 10 9 K (see curve a in figure 
2). It is interesting to point out that the existence of a “hot” equilibrium state can be 
recovered analytically by solving the energy balance equation in the limit in which Compton 
cooling dominates and using an approximated expression for the Comptonized radiation 
spectrum (see Titarchuk 1997, Titarchuk, Lapidus & Muslimov 1997). In the case in which 
Comptonization is effective and produces a rather flat spectrum with a <C 1 and the 
scattering depth at the stopping radius exceeds unity it is T ~ 10 9 K quite independently 
on the value of yo (Titarchuk, private communication). 



In the following we consider in more detail “hot” atmospheres. As we discussed, in 
this picture only a small fraction of the total luminosity (L in /<; 0.2 — 0.3) is produced 
below ijt r , where the scattering depth is already larger than ~ 10. Bearing this in mind, 
we can compute the radiation held only for y < y tr with a suitable boundary condition at 
y tr . The consistency of the computed flux at the boundary with the assumed value of L in 
will be checked a posteriori. Radiative transfer is solved by means of the characteristics 
method discussed in paper III and yields the full spatial, angular and frequency dependence 
of the photon occupation number / = c 2 / /2h 4 i/\ where / is the specific intensity. In the 
case of a static, plane-parallel atmosphere / depends on one “spatial” variable, r, on the 
cosine of the angle between the photon propagation direction and the vertical direction, y, 
and on the local photon energy E = hv. The spectrum observed at infinity is obtained 
redshifting to energies E <*, = E( 1 — R g /R *) 1 ^ 2 that at the top of the atmosphere. We 
used as boundary conditions / = 0 at the top of the atmosphere for incoming rays and a 
regularity condition (df /dr = 0) at the inner boundary for outgoing rays. The interested 
reader is referred again to paper III for all computational details about the source term 
which includes relativistic e-p, e-e bremsstrahlung and Compton scattering. In particular, 
Compton scattering is treated in its more general form by direct evaluation of the integrals 
of the Compton scattering kernel. To avoid numerical instabilities which arise when more 
accurate formulae are used (see e.g. Shestakov, Kershaw, & Prasad 1988), the Kompaneets 
approximation was retained in the calculation of the Compton energy exchange rate. 

The resulting monochromatic energy density is shown in figure 1 for L^jL e = 7x 10 -2 , 
L in = 0.2LOO and y 0 = 20 g cm -2 at different optical depths. The dashed line represents the 
spectrum observed at infinity, i.e. corrected for the gravitational redshift. This model has 
been computed solving the transfer equation for 20 values of y, 25 energies in the range 
0.017 MeV < E < 5 MeV and 10~ 3 < r < 8. The considered luminosity is close to the lower 
limit for the existence of “hot” solutions with y 0 = 20 g cm -2 , and it can be thought to be 
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well representative of situations where photon energies are large enough to produce pairs, 
since T 7 is anticorrelated with L^ (see paper I). The emergent flux, in fact, is peaked at 
~ 500 keV and is characterized by a high-energy tail. The gas and radiation temperatures, 
together with the correspondent profiles derived in paper I, are shown in figure 2. The lower 
values of T (about a factor 3) are mainly due to the different calculation of T 7 which enters 
in the energy equation. Now the radiation temperature is directly evaluated from the mean 
intensity and not obtained from a phenomenological equation, as in paper 1 . In the present 
calculation, the atmosphere appears to be nearly isothermal. In both cases the energy 
balance is practically obtained equating Compton heating and cooling, which gives T ~ T 7 . 
Now T is always close to 10 9 K and the computed luminosity at the inner boundary is 
~ 0.02 L e ~ 30%L 00 (quite close to the input value), so conduction is not expected to be 
much effective. 


3. Pair Processes 


We use the results derived in the previous section to estimate the pair number density 
n .|_, under the simplifying hypothesis that pairs thermalize with electrons and protons in the 
atmosphere. Clearly this is just an ad hoc assumption that needs to be checked a posteriori 
(see the discussion at the end of the section). 


Let us consider the process 


7 + 7 e + + e ; 


( 6 ) 


quantities relative to the two photons will be denoted with indices 1, 2. The absorption 
coefficient for a photon of dimensionless energy X\ = hui/m e c 2 propagating in direction Hi 
through a radiation held of intensity 1-2 = /(t, £ 2 ,^ 2 ) is given by 


1 r I 2 

a 77 (xi) = — — a ( x ) (1 — /a) dx2dn 2 d(J)2 , 

he J 3/9 


( 7 ) 
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where x 2 = X1X2 (1 — /i) /2 is the square of the photon energy in the center-of- 
momentum frame, /x is the cosine of the angle between the two photon directions, 0 2 
is the azimuthal angle and cr(x) is the photon-photon pair production cross-section 
(see e.g. Jauch, & Rohrlich 1976). Replacing 0 2 with 0 = 0 2 — 0i and introducing 
x\ = X 1 X 2 [1 — cos (61 ± 62 )] /2 (see Stepney, & Guilbert 1983), expression (7) becomes 


a 


77 



h 

XlX 1 


F(x + ,a;_) . 


( 8 ) 


The function 


_. . . f x + a (x) x 3 dx 

F(x + ,X-)= / 

■’ x - sj{x 2 + — a; 2 ) (x 2 — x 2 _) 


( 9 ) 


does not depend on the specific intensity so it was calculated once for all on a fixed grid of 
values of x± and stored in a matrix. The outer double integral was then computed using a 
spline interpolation. Once the absorption coefficient is known, the pair production rate _R 77 
is obtained performing two further integrations over angles and energies 


R 


77 — 



( 10 ) 


Following Svensson (1982) and introducing 6 = KT/m e c 2 , the pair annihilation rate 
can be expressed in terms of a dimensionless function A( 6 ) as 

n+ = n + ri-cr 2 e A (9) , ( 11 ) 


where 

A (9) ~ - 

1 J ~ 1 + 29 2 / In (2r)9) ’ 

77 ~ 0.56146 and r e is the classical electron radius. Combining the previous expressions, the 

pair balance equation becomes in the stationary case 



_R 77 — z{ 1 + z)n 2 cr 2 e A (9) = 0 , 


(13) 
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where z = n + /n , n is the proton number density and the charge neutrality condition 
n_ = n + + n has been used. The positive root of equation (13) is 

Z[t) = 4 (\A + n 2 A{6 ^f 2 ~ X ) ‘ (14) 

Using the specihc intensity and the run of thermodynamical variables derived in section 2, 
we calculated the proton and positron number densities shown in figure 3, together with 
z. As it can be seen, while in the inner atmospheric regions no relevant pair production is 
expected, the value of z becomes ~ 10 for r ~ 10 -3 and tends to increase in the external 
regions. 

Since the Eddington limit in a pair plasma is lowered by a factor m e /m p , there is 
the possibility that large values of z are never reached, because pairs are accelerated 
outwards by the radiative force as soon as they are created. The characteristic timescale 
for radiative acceleration t acc ~ (2 h/a) 1 / 2 ~ 5 x 10~ 9 (/r/10) 1 A(L 00 /L S ) -1 / 2 s, h ~ 10 cm 
is the scale height of the envelope and a ~ ( GM/R*)(L/Le) is the radiative acceleration, 
should be compared with the timescales for particle collisions and for the development 
of plasma instabilities. Coulomb collisions establish thermal equilibrium in a time 
t co ii ~ 6 X 10~ 7 (T/10 9 ) 3 / 2 (n/10 18 ) -1 s, while the growth time for the counter-streaming 
instability is t ins ~ 2 x 10 -11 (n/10 18 ) -1 / 2 s (see e.g. Melrose 1986). Plasma instabilities have 
therefore enough time to produce microturbulence which, in turn, will efficiently couple e± 
to the ambient plasma. Enhanced particle scattering prevents pairs from escaping and also 
allows them to thermalize with atmospheric electrons. Our initial assumption that pair 
production/annihilation is in equilibrium seems indeed justified. 
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4. Discussion and Conclusions 

The existence of “hot” solutions, first proposed in paper I, has been confirmed 
by means of a more detailed calculation of the radiation field. Our frequency- and 
angle-dependent approach allowed us to compute the specific intensity and hence to 
evaluate the photon-photon pair production. The temperature and density profiles are 
close to those presented in paper I and the spectrum peaks at about 500 keV. 

First of all, we wish to comment about the possibility of getting a “hot” solution 
started. Of course, the thermal energy stored in a “hot” atmosphere is much larger than in 
a “cold” one, and it is roughly 

U th ~ 4vr Rly tr (—) , (15) 

V m P J 

where T is the average temperature. Assuming KT ~ 100 keV, Uth ~ 10 31 erg, a factor 
~ 100 above the “cold” case. The transfer of energy from protons to electrons in our model 
is described by equation (1). Obviously, this is a rather crude estimate which does not 
enter into the details of the physical processes responsible for the energy exchange. Heating 
associated with the proton stopping might produce temperatures only in the keV range, so 
the transition between “cold” and “hot” states could require some more efficient heating 
process, like dissipation of shock waves or magnetic field reconnection. We note, however, 
that the proton bulk kinetic energy at the neutron star surface is ~ 100 MeV, which is 
much higher than the temperature of “hot” solutions. In any case, the smallness of the 
energy content of the atmosphere compared to the value of the luminosity, L <*, ~ 10 3 ' erg/s, 
strongly indicates that there should be no severe physical hindrance to drive the transition 
on time scales of microseconds. 

On the other hand, one can argue that the “hot” state should be short-lived. In fact, 
as shown is section 3, pair production is expected to become important, at least in the outer 
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layers where the pair density reaches its equilibrium value at z ~ 1 — 10. Due to the extra 
opacity produced by e ± ) the critical luminosity in the atmosphere becomes a factor 1 + 2 z 
lower than the Eddington limit. In “hot” models with L <*, ~ 0.1 Lg, where pairs appear 
to be coupled with the plasma, part of the envelope will become necessary dynamically 
unstable when z ~ 5. This conclusion probably remains valid even if the appearance of 
pairs alters the temperature profile and in a complete self-consistent solution the value of 
z is lower than that considered here. One should expect the outer layers to be expelled 
very quickly at the onset of the “hot” state: as a consequence accretion may be inhibited 
with the possible production of a relativistic shock wave. This suggestion is corroborated 
by the consideration that the momentum fluxes of the two flows are within an order of 
magnitude, as can be seen comparing M acc v acc and M out v out « M out a , where M out ~ 47r B%y 
is the mass in the unstable part of the atmosphere and a GM[{ 1 + 2 z)L/L e - 1 }/Rl is 
the acceleration. Using the numerical results of section 3 and assuming M acc ~ L(R*/GM), 
v acc = c(R g / -R*) 1 / 2 , which gives an upper limit for the momentum flux of the ingoing 
material, it turns out that M acc v acc ^ 10 M out v out . The decrease of the accretion rate and 
of the luminosity may push the system in a regime where only the “cold” solution exists, 
giving rise to an on-off behaviour. 

The “hot” state should be characterized by a spectrum close to that shown in figure 
1, with typical emission at ~ 100 keV. Moreover, pair bursts are expected, each consisting 
of ~ 47 T Rlyz/rn p ~ 6 x 10 35 particles; an upper limit for the integrated luminosity in the 
annihilation line will be of order of 10 29 erg. These hot, optically thin and Comptonized 
envelopes bear some resemblance with hot coronae above accretion disks, although the 
geometry, the heating process, and the triggering of the pair instability may be different. 
The onset of the “hot” state represents a possible physical mechanism for producing 
efficiently high energy radiation from weakly magnetized, accreting neutron stars and may 
be of interest in connection with hard X-ray transients, although we are aware that some 
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simplifying assumptions, like spherical symmetry and absence of magnetic fields, may limit 
the applicability of our results. 

At present hard emission (E ^ 30 keV) has been observed from a number of sources, 
believed to contain NSs, with SIGMA and BATSE (see e.g. Tavani, & Liang 1996; Vargas 
et al. 1996). At least seven X-ray bursters, including Aql X-l, and other LMXBs show 
transient emission in the range 30-200 keV, with intensity anticorrelated with that in the 
soft band. These observations are of the utmost importance, since hard X-rays emission 
were previously associated only with X-ray pulsators and black hole candidates. Tavani, & 
Liang have proposed a possible explanation for the transient, hard state of these sources 
in terms of magnetic reconnection in the inner regions of an accretion disc around the NS, 
but the nature of the primary emission mechanism is still an open issue. The total X-ray 
luminosity of these sources is ~ 10 36 — 10 37 erg/s, close to that considered in our model, 
although the duration of the hard state is several days. In our picture, this would imply 
that the appearance of e 1 * 1 " pairs, even if it affects the accretion process, does not induce the 
transition to the “cold” state. Hard states of shorter duration would be impossible to detect 
with BATSE and SIGMA, but they may be observed by the instrumentation on board of 
XTE or SAX. 

We are indebted both to an anonymous referee for drawing our attention on the 
importance of thermal conduction in “hot” models and to the second referee, Lev Titarchuk, 
for his clarifying comments. The revised version of this paper greatly benefited from 
discussions with Fred Lamb, Luciano Nobili and Luca Zampieri. We are also grateful to 
Marco Tavani for useful discussions about hard X-ray transients. 
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Fig. 1.— Monochromatic mean intensity for the model Loo = 7 x 10 2 L e , ?/o = 20 g cm 2 ; 
different lines correspond to equally spaced values of logr in the interval [—3,0.9]. The 
emerging redshifted spectrum is also shown (dashed line). 

Fig. 2.— a) The gas temperature (full line) and the radiation temperature (dashed line) for 
the model of figure 1. b) Same results from the frequency-integrated analysis of paper 1. 

Fig. 3.— Proton (full line), positron (dashed line) number densities (in units of 10 22 cm -3 ) 
and z = n + /n (dash-dotted line) versus the scattering optical depth for the model of figure 
1 . 
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